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CONTROL OF SPATIALLY HETEROGENEOUS AND 
TIME-VARYING CELLULAR REACTION NETWORKS: A NEW 

SUMMATION LAW 

MARK A. PELETIER, HANS V. WESTERHOFF, AND BORIS N. KHOLODENKO 

Abstract. A hallmark of a plethora of intracellular signaling pathways is the 
spatial separation of activation and deactivation processes that potentially re¬ 
sults in precipitous gradients of activated proteins. The classical Metabolic 
Control Analysis (MCA), which quantifies the influence of an individual pro¬ 
cess on a system variable as the control coefficient, cannot be applied to spa¬ 
tially separated protein networks. The present paper unravels the principles 
that govern the control over the fluxes and intermediate concentrations in spa¬ 
tially heterogeneous reaction networks. Our main results are two types of the 
control summation theorems. The first type is a non-trivial generalization of 
the classical theorems to systems with spatially and temporally varying con¬ 
centrations. In this generalization, the process of diffusion, which enters as the 
result of spatial concentration gradients, plays a role similar to other processes 
such as chemical reactions and membrane transport. The second summation 
theorem is completely novel. It states that the control by the membrane trans¬ 
port, the diffusion control coefficient multiplied by two, and a newly introduced 
control coefficient associated with changes in the spatial size of a system (e.g., 
cell), all add up to one and zero for the control over flux and concentration. 

Using a simple example of a kinase/phosphatase system in a spherical cell, we 
speculate that unless active mechanisms of intracellular transport are involved, 
the threshold cell size is limited by the diffusion control, when it is beginning 
to exceed the spatial control coefficient significantly. 


1. Introduction 

Extracellular information received by plasma membrane receptors is processed, 
encoded, and transferred to the nucleus through activation and spatial relocation of 
multiple signalling components. Receptor activation triggers signalling responses 
associated with the mobilization of a plethora of adapter and target proteins to 
the plasma membrane (Haugh and Lauffenburger 1997; Kholodenko et al. 2000b). 
Proteins activated e.g. by phosphorylation at the cell surface travel to stimulate 
critical regulatory targets at various cellular sites including the nucleus. During 
and after the travel, these proteins are inactivated, e.g. through dephosphorylation. 
For instance, a protein phosphorylated by a plasma-membrane associated protein 
kinase can be dephosphorylated by a phosphatase in the cytosol or nucleus. The 
transport between cellular locations, where the activation and inactivation occur, 
is passive (thermal diffusion), driven by spatial gradients in the concentrations. 
In some cases, not to be discussed here, motor proteins or endocytosis may be 
involved (Kholodenko 2002). Previous work estimated that these spatial gradients 
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may be large, and that therefore diffusion may contribute to the control of signal 
transduction (Brown and Kholodenko 1999; Kholodenko et al. 2000a). 

For spatially homogeneous reaction networks, the control over fluxes and in¬ 
termediate concentrations has been studied both experimentally and theoretically 
(reviewed in (Fell 1997)). The control is generally quantified as the extent to which 
any type of molecular process influences a system variable, such as the flux or 
concentration. The control coefficients are defined as the ratios of the fractional 
changes in the system variable to that of the biochemical activity which caused the 
system change. Mathematically, the stability of a system steady state is required 
and the changes are considered as infinitesimally small, i.e., the coefficients are ex¬ 
pressed as the log to log derivatives (Fell 1997). Important principles underlying the 
control of biochemical reaction networks in well-stirred reactors have been worked 
out recently. For instance, it can be shown that a large increase in the activity of 
a single enzyme does not result in a substantial increase in the flux for almost any 
metabolic network; marked flux increase can be achieved by a concerted modula¬ 
tion of several pathway reactions (Kacser and Acerenza 1993). In metabolic control 
analysis (MCA), this result is related to the so-called summation theorem, which 
states that the sum of the enzyme control coefficients adds up to 1 (Kacser and 
Burns 1973; Heinrich and Rapoport 1973; Westerhoff and van Dam 1987). There 
is no rate-limiting enzyme, the control is shared between all network processes. 

The present paper unravels the principles that govern the control pattern in spa¬ 
tially heterogeneous cellular networks. We demonstrate new properties of control 
that result from the spatial aspect of diffusion and separation of signalling reactions. 
The control summation theorems, relevant for these networks, are formulated and 
proven. 


2. Methods 

2.1. Description of the system. We consider living metabolic systems consisting 
of n internal chemical species Yi,... ,Yn, with concentrations ci,...,c„. These 
species interact via m reactions which can be represented as 

n n 

^ E^b^- j = (1) 

i=l 

where are the dimensionless forward and reverse stoichiometric coefficients 
for the consumption of species Yi in the forward and reverse reactions j. We 
write Nij = for the net stoichiometric component. Reaction j has 

rate avjVj{ci, ..., c„), so that in the absence of diffusion the time evolution of the 
concentrations c, should be given by 

d 

— c, = ^Nijay.Vj, i = 1,... ,n. 

Ci has the dimension of concentration [M], or mole number per unit volume. The 
rate of the reaction 1, i.e. Vj, has the dimension of concentration per unit time [M/s]. 
It is a function of the kinetic parameters of the enzyme catalyzing that reaction (or 
if non-enzymatic of the chemical reaction itself), e.g. Km, the concentration of the 
enzyme, and the concentrations of substrates, products and allosteric modifiers of 
the reaction. The dimensionless parameter is the tool we shall use to modulate 
the activity of reaction vj. 
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In this paper we will be interested in cases in which concentrations can not be 
assumed homogeneous in space inside the systems. For thermodynamic reasons, 
living systems need to be open for at least some chemical compounds. These sub¬ 
stances are transported across the system membrane, usually by transport proteins, 
i.e. effectively there is a sink or source for some of the species on the system bound¬ 
ary, e.g. the (plasma) membrane. The spatial separation between the sink or source 
and the chemistry that takes place in the bulk phase then creates spatial variations 
in concentration. 

When the system is not homogeneous and not at steady state, the concentrations 
need to be specified as functions of time and space (i.e. Ci{x, t)), and boundary fluxes 
are required. The principal difficulty accompanying the spatial inhomogeneity of the 
concentrations of any component Yi, which makes the classical MCA inapplicable, 
is that the number of ‘explicit’ variables associated with Yi becomes infinitely large. 
The reactions are assumed to take place in an enclosed three-dimensional space that 
is represented by a bounded set of xi, a ;2 5 a^ 3 -coordinates, i.e. the volume F C 
For instance, V can correspond to the cytoplasm or the mitochondrial matrix. 
Inside V the species react with each other, according to Eq. 1; in addition they 
undergo diffusive transport. The concentrations c* of the species are functions of 
time (0 < t < tend) and space (x G F), and their evolution is governed by the 
balance (reaction-diffusion) equations (Katchalsky and Curran 1965) 

d 

^ ^ ^ X G F, 0 < t < tend- (2) 

The parameter Di is the diffusion coefficient of species F. It is considered homo¬ 
geneous. For Cartesian coordinates the Laplacian A is defined as 

^ ^ ^ ^ 

“ dxj ^ dxl ^ dxj ■ 

The operator A has the dimension of 1/length^ Di has the dimension of 

surface per time [m^/s] and is assumed to be the same throughout volume F. In 
order to assess the importance of diffusion of a substance for the behavior of the 
system, we shall modulate the diffusion processes. To this aim, Di is multiplied 
by the dimensionless modulation parameter ajj., which has the value of 1 in the 
reference state. 

In addition to Eq. 2, which holds in the bulk, represented by F, we need to 
specify the behaviour at the boundary S of the set F. Eor a species F for which 
the membrane is a passive barrier we assume a closed boundary condition, 

dc ■ 

-ajj.Di-^ = 0 on A. (3) 

The vector u is the normal vector of unit length pointing out of F, and dcijdv is 
the derivative of Ci in the direction of v. This makes the quantity —aj^.Didcijdv 
the outward diffusive flux of species F- 

If a species is transported across the boundary of F, then the boundary condition 
takes a similar form, but with non-zero right-hand side, 

Qc' 

~°^DiDi-^ = afj, on S, (4) 

where fi is the transport rate per unit surface area, fi can be due to active transport 
or a (bio)chemical reaction at the membrane surface and may be a function of the 
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concentrations ci,..., Cn. There is at least one non-zero flux a/^/i (one of which 
we conveniently assume to correspond to i = 1). Like other a’s, aj^ is a parameter 
that we shall use to modulate the transport activity. This modulation is equivalent 
to changes in the activity or concentration of membrane enzymes or carriers. We 
shall refer to the ‘flux J through the system’ as the magnitude of the net export of 
Yi, integrated over the surface of V and normalized by the surface area: 

We normalize the total flow with respect to the corresponding surface area to obtain 
a surface-averaged flux. The difference in role of a reaction rate in the bulk phase 
(vj) and a transport rate (/i) is underlined by the difference in the dimensionality— 
the former is given in volume concentration per time unit (e.g., M/s), whereas the 
latter is expressed in surface concentration (density) per time unit (e.g., M • m/s). 
We comprise Eq. 3 into Eq. 4 by allowing fi to be zero. Some reactions that involve 
cytoplasmic substances only occur in the membrane or at the membrane surface. 
This is equivalent to the substrate of the reaction being exported and the product 
of the reaction being imported, and will be treated as such. 

2.2. Parameter Modulation and Control Analysis. We have explicitly intro¬ 
duced dimensionless modulation parameters through which the 

system can be modified. For the formulation and interpretation of control theorems 
associated with diffusion we need to introduce two additional forms of modulation, 
associated with space and time, in the following way. We assume that the spa¬ 
tial and temporal variables Xi,X 2 ,X 3 , and t, are related to a reference spatial and 
temporal frame Ci)C 2 )^ 3 j and via dimensionless modulation parameters az, and 
at'. 

Xi = aL^i, X 2 = aL^ 2 , X 3 = aL^ 3 , and t = atT. (6) 

Note that a reference time interval [0, Tend] is rescaled to an actual time interval 
[0, OtTend] = Tend]; and similarly, the reference volume V is scaled to an actual 
volume olV (Figure 1). 



Figure 1. The domain V is scaled by geometric multiplication 

We shall be explicit about the cause-and-effect relationship that exists between 
parameters and state. The parameters are the kinetic constants in rate functions Vj 
and fi and diffusion coefficients Di, together with their modulators aj^., and 
af^. We also consider the volume V and the time interval [0,Tend] as parameters, 
together with the space and time modulators ar and at- A final ‘parameter’ that 
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we need to take into account is the distribution of all species in space at the initial 
time. 

For a given choice of all parameters the concentration of each of the species Yi 
is determined uniquely as a function of space and time, according to Eqs. 2 and 4. 
The choice of parameters yields a a unique set of functions Ci{x, t), each defined at 
every point x G ulV and at each time t G [0, OtTend]- 

Note that the functions Ci are defined on different domains for different ol and 
at ■ For comparison purposes it will be more convenient to think of Ci as functions of 
the reference coordinates ^ and r and the modulation parameters a, i.e. Ci(^, r, a), 
where ^ GV, t G [0,Tend], and a = az,, at). The functions Ci{^,T,a) 

determine the state of the system associated with the choice of parameters. 

The change of coordinates, from {x,t) to induces a small modification in 

Eqs. 2 and 4, which now become 


d "" 

O^t ^ ^ 1 ^ ^ 0 < T ^ Tend (7) 

^ i=i 

dc- 

- a^^cxjy.Di-^ = ajji ? G S', 0 < t < Tend (8) 

Note that the differentiation in the operators A and d/di^ is now with respect 
to f. The definition of the flux, however, remains unchanged under this change of 
coordinates, because of the surface averaging: 


In the vein of Metabolic Control Analysis we will consider a reference state 
associated with the choice = a^ = a^ = 1 and a second state, 

resulting from an infinitesimal change in one of these modulators. The difference 
between the two states will be taken to characterize the control that the modulated 
parameter has on the state of the system. The dimensionless parameter is 
the tool we shall use to modulate the activity of reaction Vj. Departures of 
from 1 correspond to a modulation of enzyme catalytic activity or concentration, 
or to a proportional modulation of all rate constants of the reaction step if it is not 
enzyme catalyzed (this modulation leaves the equilibrium constant unchanged, for 
details see (Kholodenko and Westerhoff 1993, 1994, 1995), where this approach was 
introduced). Eor instance, if reaction 1 is a protein-protein interaction (a central 
reaction in signaling networks), the rate Vj depends on two kinetic constants, fcon 
and fcoff, and a change in a„^. corresponds to equal proportional change of both 
constants that leaves the dissociation (equilibrium) constant unchanged. The 
ratio of the fractional change in a state variable of the reaction network, such as 
J or Cfc, and the fractional change in a„^. determines the control coefficient with 
respect to reaction rate Vj (in the limit of infinitesimally small changes): 


ai = 


din J 


din av 


= 


din Ck 


din av 


where ay^ = 1 corresponds to the reference state (all other parameters are assumed 
fixed). Similarly, considering two states corresponding to small modulations in di¬ 
mensionless parameters, a/^ and a^., we define the control coefficients with respect 
to transport reaction and diffusion, Cj. and . 
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An intrinsically novel control coefficient emerges as the heterogeneous spatial 
organization of cellular reaction network is taken into account. This coefficient is 
obtained when the size (L) of the system is modified through a modulation of the 
dimensionless parameter as shown schematically in Fig. 1, whilst keeping hxed 
all other parameters. In a cell setting, these parameters include the volume and sur¬ 
face concentrations of enzymes and other molecular forms, the total concentrations 
of which are conserved (for instance, the total amount of phosphorylated and un- 
phosphorylated forms of a protein) in network reactions described by Eqs. 7 and 8. 
Therefore, increasing the cell size also implies adding conserved chemical moieties 
and additional enzymes both to the bulk aqueous phase and cell membranes. Ac¬ 
cordingly, the ‘spatial’ control coefficient is defined as Cj^ = din J/c?Inwhere 
only the spatial variables are modulated, as described above (Eq. 6). It is of note, 
that the parameter L can be interpreted as the characteristic cell size, e.g., the vol¬ 
ume to surface ratio, whose fractional changes are equal to the fractional changes 
in a^. 

A special case of Eq. 7 arises when we only consider steady-state, non-equilibrium 
configurations; these correspond to solutions ci(^),..., c„(^) of the steady-state 
equations 

m 

-al'^ajy.DiAci = G E, (10) 

i=i 

with boundary conditions 8. For the purposes of this paper we assume that solutions 
of this set of equations are locally unique and asymptotically stable. 


3. Results: the first summation theorem 


In a metabolic network where concentrations are homogeneous, the sum of the 
flux control by all biochemical reactions equals 1 (Kacser and Burns 1973; Heinrich 
and Rapoport 1973). We generalize this summation theorem for systems where 
concentrations may not be homogeneous and where transport occurs — here the 
sum over all the reactive processes need not be 1. The flux may be partially 
controlled by diffusion or by transport and it is because of this control that the 
classical summation theorem no longer holds. The generalized summation theorem 
takes the form 

n m n 

+ = ( 11 ) 
i=i j=i i=i 

In words this theorem reads: the total control by all diffusion, reaction, and trans¬ 
port processes on any steady state flux equals 1. The proof of Eq. 11 is based on 
Euler’s theorem on homogeneous function, which can be stated as follows. Let g 
be a function of pi,... ,p„, such that for all pi and for all A > 0, 

5(A^ipi,..., X^’^pn) = X'*g{pi, ...,p„), (12) 


for some /3i,... G K., i.e., g{pi, ■ ■ ■ ,Pn) is the same function after the trans¬ 

formation 


g = X^g, Pi = A^'p*. 


(13) 


Then 


^ 91n|g| 


= 7- 


(14) 
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It may be noted that for some of the parameters pi, Pi may equal zero. These 
parameters are not modulated in Eq. 13, and Eq. 14 shows that these parameters 
are absent from the summation theorem. 

For 7 and all Pi equal to 1, Eq. 14 reminds of the flux control summation theorem 
of Metabolic Control Analysis (Kacser and Burns 1973; Heinrich and Rapoport 
1973), but it may not be clear a priori what function and which parameters should 
be considered. When searching for other than the traditional summation theorems, 
a strategy may prove useful. One strategy is that of leaving the system in essence 
in the same state (cf. (Kacser and Burns 1973)). In terms of the system under 
consideration this translates into the concentrations Ci (^, r) remaining the same 
when the parameters pi of Eq. 13 are changed. 

In order to implement the above strategy and examine under what type of pa¬ 
rameter changes Ci(^, r) might be constant, one may inspect the equations that 
define c^. For steady states Eqs. 10 and 8 define cpp). Multiplying both equa¬ 
tions by the same factor A should leave their solutions unchanged. Consequently, 
multiplying all and all af^ by that same factor should not change 

the solutions cpp) either: the terms in the equations are homogeneous functions of 
first order of the parameters and af^. When only these parameters are 

modulated, the specific transport rate remains unaltered, because the metabolite 
concentrations remain constant and neither q;„^. , or occurs in the function 
/i(ci,..., c„). Consequently fi remains unchanged, and therefore we can calculate 
the corresponding flux as a function of A, 

•/(A) = A^ /^/i = AJ(l). 

We can reformulate this result as follows: we have shown that the flux J, as a 
function of the parameters , and a/^, satisfies 

J (Aaj,..., Xa Dn ) Aat,j,..., Xcxy ^, \a ,..., Xa ) 

= XJ{aDi , ■. •, Q;_d„, ..., ,..., «/„). 

Equality 11 now follows from Eq. 14. 

It should be remembered that the parameters , ay^ , and a represent the ac¬ 
tivities of diffusion, chemical, and transport processes, respectively. Consequently, 
Eq. 11 really refers to a theorem concerning the control by diffusion, chemistry and 
transport. It is an extension, both in concept and in proof, of the familiar flux con¬ 
trol summation theorem of Metabolic Control Analysis. The latter only dealt with 
the middle terms of the left-hand side and then only the enzyme catalyzed reac¬ 
tions thereof, it assumed concentrations to be homogeneous in space, and neglected 
explicit transport terms. 

A further extension to the above methodology will generalize the above sum¬ 
mation theorem to systems that depend on time: Dropping the assumption of 
steady state, Eqs. 7 and 8 become the equations that define Cj(^,r). Multiplying 
both these equations by the same factor A again leaves their solutions unchanged. 
Consequently: 


-C, 




i=i 


Eg. 


J(t) 


= 1 . 


(15) 



8 


MARK A. PELETIER, HANS V. WESTERHOFF, AND BORIS N. KHOLODENKO 


Here coincides with the time-control coefficient defined in (Acerenza et al. 

1989). How can we understand this theorem and in particular this new time-control 
coefficient? Let us first consider the situation that would ensue from simultaneously 
and equally increasing all the parameters ay^, and by 1 % in the absence 
of an increase in at- Because all activities increase, also the flux J should increase. 
However, the system should also change more quickly in time, i.e. the increase in 
magnitude of J should occur earlier. To obtain a proportionate increase in J, one 
should look at an earlier point in time, i.e. time should be 1 % earlier, or, keeping 
in mind that t = atT, at = 0.99. 

should not be confused with the time-dependent control coefficients de¬ 
fined by (Westerhoff and van Dam 1987; Heinrich and Reder 1991). For a further 
discussion of the distinctions, see also (Kholodenko et al. 1997) for the control 
analysis of relaxations in the vicinity of the steady state. 

Often, as time goes to infinity, the system relaxes to a steady state. The control 
by time on the flux then reduces to zero and the summation theorem for the time 
dependent control coefficients reduces to the one for steady state control coefficients 
(Eq. 11). 


The summation theorem 15 is more general than suggested by the above deriva¬ 
tion. It can be derived without requiring explicit equations for the time evolution 
of the system, such as Eqs. 7 and 8. The proof then considers an actual physical 
transformation of the system under consideration: all the time dependent elemental 
processes are increased by the same factor A. This implies that everything happens 
in the same way, but faster by that same factor A. We represent this concept of 
‘faster’ by the parameter at, as above (since t = atT, at < 1 corresponds to ‘faster’, 
and at > 1 to ‘slower’). Writing the parameter-dependence explicitly. 


we have 






or in other words, c* is homogeneous of zeroth order in a^., a„^, a/^, and 1/at. 
Using Eq. 14 this yields the summation theorem for concentrations. 


n m n 

-cr + ^ c/,/ + ^ ^ c;/ = 0. (16) 

i—1 j—^ ^=1 

where Ck is short for Ck{^,T). A similar result for the flux, 

J , Aa , Xayj, \aj^, — \J{T,a , a^^ , a at), 

shows that the flux is homogeneous of first order in a^., a^^, O'f^, and 1/at, resulting 
in the summation theorem Eq. 15. 


4. A NOVEL SUMMATION THEOREM 

The above theorem is an extension of the theorem of (Acerenza et al. 1989) to 
cases with transport and diffusion. The theorem we shall develop now is intrinsi¬ 
cally new as its origin lies in spatial aspects that have hitherto been neglected in 
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Metabolic Control Analysis (Fell 1997). The strategy we follow is again to con¬ 
sider a transformation for which the system properties should not change. We first 
discuss the nature of diffusion, as illustrated by Einstein’s diffusion equation: 

6x = a/ nijDdt, 

where nrr is a number of order 1 that depends on the dimensionality. In words, 
the displacement due to diffusion increases with the square root of the diffusion 
coefficient. Consequently, increasing all diffusion coefficients by should increase 
all displacement by the factor A. 

This becomes obvious when considering Eq. 7. Replacing a^. by A^a^. and 
ctL by Atti leaves Eq. 7 unchanged. The boundary condition 8 in a natural way 
suggests replacing af^ by Aa/^, so that this equation is also essentially unchanged. 
As a result this choice of parameters leaves the solution Ci invariant: 

(C? ^ Di 5 ^ 5 (^5 U, O; , Oivj : j CX L ^ CXt'). (1^) 

We immediately deduce the second summation theorem for concentrations, 

n n 

2E^^z+EC'/z+C'£ = 0, (18) 

i—1 i—1 

where c stands for any Cj (^, r). Similarly, since the flux satisfies 

J {t ^ \ CX , ^Vj : : Aoi^ lO^i) — \ J {t ^ (X , CXvj : ^ fi : ^ L : ^ 

the associated summation theorem results 

n n 

2 E + E 

i=l i=l 

Naturally, if we only consider steady states, then Ci{^, t) reduces to Ci(^), and J(t) 

to J. 

The modulation of V, via implicated in the evaluation of the corresponding 

control coefficient is one in which the size of the system changes; it is important, 
however, that the system size changes whilst keeping the volume concentration of 
enzymes in the bulk aqueous phase constant and the surface concentrations of the 
membrane enzymes constant. 


To illustrate the modulation introduced above, consider the case of a one-dimensional 
spatial domain of the form (0, oo); the membrane is thought to be at a; = 0, and 
the cell is deemed so large that the ‘other end’ is effectively at infinity. For such a 
‘half-infinite’ cell, in non-equilibrium steady-state, Eq. 19 reduces to 

n n 

2EC'D.+EC'h = 1 (20) 

where J is the flux at the border {x = 0); C'l vanishes as a result of the infinite 
size. Figure 2 shows how the solution Ci(x) depends on the parameter modulation 
by the factor A given in Eq. 17. Note that here Ci is plotted as a function of x; 
when plotted as a function of the graphs for different A coincide (cf. 17). 
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Figure 2. The scaling in A for a two-species system 


phosphatase 



^ _kinase 





Figure 3. Geometry 1: a slice bounded by parallel membranes 



kinase 
at the 
membrane 


Figure 4. Geometry 2: a spherical cell 
5. An explicit example system 

For certain simple cases control coefficients can be calculated explicitly. We do 
so here for a kinase/phosphatase example that we also discussed in (Brown and 
Kholodenko 1999; Kholodenko et al. 2000a). There are two species, Y and YP, 
which are the unphosphorylated and phosphorylated form of a given protein; the 
kinase converts Y into YP, and the phosphatase does the reverse. 

Besides providing an illustration of the concepts and results discussed above we 
wish to demonstrate the influence of the system geometry on the control coefficients. 
We consider two cases: 

(1) V is bounded by two parallel membranes; the kinase reaction is localized 
to one membrane, the phosphatase to the other (Figure 3); 

(2) F is a spherical region (e.g. a cell or an organelle) bounded by a membrane; 
the kinase is localized at the membrane, but the phosphatase is distributed 
throughout the bulk (Figure 4). 
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As we shall see below, in the case of parallel membranes the size control coefficient 
C/ is negative: an increase in the distance between the membranes reduces the flux 
through the system. However, when the phosphatase is not membrane-bound but 
distributed, as in the second case, increasing the system size actually increases the 
flux: C'l is positive. 


5.1. Two parallel membranes. We assume that the system is large with respect 
to the distance L between the membranes and therefore adopt a one-dimensional 
formulation. We also restrict the study to (non-equilibrium) steady state. 

Since no reaction takes place in the bulk, the equations satisfied by Y and YP 
in the bulk (Eq. 10) reduce to 

-al^ajjD YP" = -al^ajjD T" = 0, 0 < C < A- (21) 

Note that the spatial variable is as discussed above; the physical variable x is 
given by a; = 

For the kinase and phosphatase reactions we assume the rate functions 


akkk{KkY — YP) and apkp(YP — KpY). 


Note the dimensionality of kk,p: the rate functions are surface fluxes, but YP and 
Y are bulk concentrations; kk and kp therefore have dimension length/time. The 
coefficients ak,p and Kk,p are dimensionless. 

The kinase and phosphatase reactions both enter the description as boundary 
conditions (cf. 8): 


ct^^ajjDYP' = apkpiYP — KpY) 

a^^ctoDY' = —apkp{YP — KpY) 


at ^ = 0 (phosphatase) (22) 


-al^ajjD YP' = -akkk^KkY -YP) 

, a.t £ = L (kinase) 

-al^ajjD Y' = akkk{KkY -YP) 


(23) 


Equation 21 implies that the spatial gradients of YP and Y (i.e. YP' and Y') 
are constant in £; let us set 6 = YP'. From the solution of equations (21-23) the 
relevant information for our purposes is 6: 

^ ^ M{f3pakKk - (3kapKp) 

+ /3p) + L(3k(3p 

Here M is the total concentration Y + YP, which is independent of £ by Eq. 21, 
and 

f3k = ctkkkil + Kk) and (3p = apkp{l + Kp). 

Because of the spatial separation of the reactions the flux as defined by Eq. 5 
equals the diffusive flux, 

J = -al^aoD YP' = -al^ajjD S. (25) 


From Eqs. 24 and 25 we can determine the control coefficients of the flux with 
respect to the kinase and phosphatase reactions and to diffusion, as well as the 
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‘length’-control coefficient:^ 

Df3p 

D{Pk + Pp) + LPkPp 
PPk 

D{Pk + Pp) -\- LPkPp 

_ LPkPp _ 

D{Pk + Pp) + LPkPp ’ 

_ LPkPp _ 

D{Pk + Pp) + LPkPp 

Clearly these expressions satisfy the summation theorems Eqs. 11 and 19: 

Cd + L-k + C'p = 1, 

2Ci + Ci + C^ + Ci = l. 

Note that both the kinase and the phosphatase reaction are boundary effects, and 
therefore give rise to boundary control coefficients of Eqs. 11 and 19. 

For concentration control we choose to consider not a pointwise concentration, 
but the difference between concentrations at opposite ends: 

YP{L) - YP{0) = LS. 


ci = 

ci = 
ci = 


Denoting the control coefficients with respect to this quantity as C"^, we have 

Dpp 


_ ___ 

D{Pk + Pp) + LPkPp 

QC _ _ PPk _ 

^ P{Pk + /dp) + LPkPp 
P{Pk + Pp) 


ch = 


P{Pk + Pp) + LPkPp 
fjc _ P{Pk + Pp) 

^ P{Pk + Pp) + LPkPp 

An example of the values of these control coefficients 




5.2. A spherical cell. For the second example we consider the same system as 
above, but now in a spherical cell of radius L. We assume spherical symmetry 
throughout. 

The kinase reaction is again localized at the membrane, and therefore results in 
a boundary condition 

-al^ajjPYP' = akkk(YP-KkY), 


-a^ ajyPY' = -akkk{YP - KkY), 


at 5 = L (kinase). 


(26) 


The spatial variable is again but it now represents a radial coordinate. The 
phosphatase reaction, which in this example is distributed throughout the bulk, 
now enters as a right-hand side in the differential equation for Y and YP (cf. 10): 

-afa^Pr^ieYPy = -apkpiYP - KpY), ^ ^ 

-afaaPrHfy'Y = apkp{YP-KpY), ^ ^ 


^Note that these formulas contain no a’s, since they result from evaluation at = 

ak = ap = 1; consequently, Pk,p = + i^k,p)- 
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L {fJ,m) L {pLm) 

Figure 5. Values of control coefficients as a function of L for the 
first example (parallel membranes at distance L, the kinase and 
phosphatase are localized at different membranes; D = 1 /s, 

kp = kk = 1 /s, Kp = 0 . 1 , Kk = 10 ). 


Again note the dimensions of kk and kp'. similar to the previous example, kk has 
dimension length/time; but here kp converts a bulk concentration into a bulk flux, 
and therefore kp has the (more familiar) dimension 1/time. As before, the total 
concentration is independent of and denoted by M: 

M = YP + Y (27) 

As a first step we investigate different boundary conditions, where we prescribe 
the flux J: 

-al^aoDYP' = J = al^ajjDY', (28) 

at ^ = L. For this choice of boundary conditions we can solve the system exactly. 
Set u = YP — KpY ; then u satisfies the equation 

= —apkp{l + Kp)u, 

with boundary condition 

—ay^ajyD u'{L) = (1 + Kp)J. 

We find that u is given by 


u{0 = -J 


1 A Kp L‘^ 

ay^ajjD yLcoshyL - sinhyL 


sinh 7 ^ 


(29) 


where the inverse length scale 7 is defined by 
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It follows that 


ym = Y^M+ 

1 + Kn 


1 




+ Rp 

1 


■M - 


1 + Kp 

1 


u{C) 

m(C). 


( 30 ) 


1 + Kp 1 + Kp 

To return to the boundary conditions 26 we equate the flux J in Eq. 28 with the 
right-hand side of Eq. 26: 


J = akkk(YP{L) - KkY (L)) = akkk ( 7 - 7 -^ 


1 -t- Kjj 


Kp Kk 
1 -t- Krn 


M 


using the expression 29 for u this yields 

d = —-7, v = akkkM— -, 

1 -I- P I + Kp 


(3 = akkk 


ajjD yLcothyL-l 


(31) 

(cotha; is the hyperbolic cotangent, cosha;/sinha:). 

As in the previous example we thus obtain explicit formulas for the control of 
the various parameters on the flux: 


ci 

ci 

ci 

ci 


1 

1 + /? 

1 (3 yLcothyL — 2-f 7 ^L^(coth^ 7 T — 1) 

2 1-\-f3 7 Lcoth 7 L—1 

1 j3 7 Lcoth 7 L — 7 ^L^(coth^ 7 L — 1) 
2\-\-f3 7 Lcoth 7 L—1 

f3 1 — 7 ^L^(coth^ 7 L — 1) 

1+ (3 7 L coth 7 L — 1 


As an example of the concentration control values we choose a point at distance 
L/2 from the center of the cell. Denoting the control coefficients with respect to 
this quantity as we calculate 

Cl = t^ci 

Cd = P ^C'ia ~ 

c; = P (c^ - 

Cl = tx{Ci-<i,+ l) 


Here 

u{L/2) 

^ KpM + u{Lf2) 


and 


2 r2 


YL 


7 L cothyL — 1 



An example of the values of these control coefficients is given in Figure 6 . 


6. Discussion 

Until fairly recently, biochemistry focused more on the time dependent aspects 
of processes than on the spatial aspects. Exceptions were the membrane-mediated 
subcompartmentation of cellular metabolism. The cellular subcompartments were 
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Figure 6. Values of control coefficients as a function of L for the 
second example (a spherical cell with a kinase on the membrane 
and a phosphatase in the cytosol; D = 1 kp = kk = 1/s, 
Kp = 0.1, Rfc = 10). denotes control of the concentration mea¬ 
sured halfway between the center and the membrane. 


mostly considered to be homogeneous in terms of concentrations and the hetero¬ 
geneity between them was analyzed in terms of the activity of transport catalysts in 
the membranes (Westerhoff and van Dam 1987). Inhomogeneity of metabolites and 
even ions within aqueous subcompartments has often been proposed (Kell 1979) but 
calculations have shown that at least for central metabolic routes aqueous diffusion 
should be fast on the time scale of the catalytic turnover of the enzymes also given 
the proximity of the enzyme molecules due to the small sizes of cells (Westerhoff 
and van Dam 1987; Westerhoff and Welch 1992). Metabolic Control Analysis was 
developed from this perspective and proved a useful way to rationalize the study 
of regulation and control of intermediary metabolism (Fell 1997) and free-energy 
transduction (Westerhoff and van Dam 1987). More than an analysis procedure. 
Metabolic Control Analysis also provided biochemists with a set of laws (theorems) 
that govern control of fluxes and concentrations in metabolic networks. One type 
of these were the summation theorems, stating that the sum of the control by all 
the individual enzyme catalyzed reactions on flux and concentration (Kacser and 
Burns 1973; Heinrich and Rapoport 1973), oscillation amplitude (Kholodenko et al. 
1997) or relaxation time constant (Heinrich and Reder 1991) should equal 1, 0, 1, 
and 1 respectively. 

Those were the days where biochemistry focused on metabolism. However, only 
a minor fraction of the known genomes encodes primary metabolism. Much coding 
capacity is devoted to regulation, in part through signal transduction (Kell 1979; 
Kholodenko et al. 2000c). Signal transduction has (i) the spatial aspect that many 
signals arrive at the cell’s plasma membrane and have to be transferred to its 
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nuclear DNA, (ii) the added complication that much of the signal is transferred 
by proteins in direct phosphoryl transfers, rather than with small molecules as 
messengers, and (hi) the feature that those proteins are often present at a thousand 
times lower concentrations than the usual metabolite. All these three features work 
in the direction that might make diffusion of the components limiting in signal 
transfer. Indeed, it was calculated that for realistic parameter values, concentration 
gradients in signal (dehned as in (Francke et al.)) should arise (Fell 1980; Brown and 
Kholodenko 1999; Kholodenko et al. 2000a), which might even force the cell to take 
refuge to alternative mechanisms of movement of the signal proteins (Kholodenko 
2002 ). 

With the realization of the possible importance of transport and diffusion for 
cell biochemistry, there came a need to update Metabolic Control Analysis so as 
to include the aspects of diffusion. This is what the first part of this paper accom¬ 
plished for the summation theorem. The result was rather simple, to the sum of the 
enzyme control coefficients also coefficients for the control by the diffusion of the 
various species needed to be added, as well as control coefficients for the transport 
processes and for the non-enzyme-catalyzed chemical reactions, before it amounted 
to 1 for flux control and 0 for concentration control. This result is perhaps not 
surprising, but has been overlooked effectively in the literature for quite a while. In 
studies of the control of growth rate of E. coli by the lactose permease, it was found 
that a sum of 1 could not be found and it was suspected that part of the control 
resided in diffusion across the outer membrane (Jensen et al. 1993; Dykhuizen et al. 
1987). It was found that none of the glycolytic enzymes had significant control on 
the glycolytic flux in yeast. This was considered a problem until strong indica¬ 
tions arose of substantial control in glucose transport into the cells (Reijenga et al. 
2001). Thus, if anything, our extension of the summation theorem to include the 
transport steps may open eyes to the possibility that transport controls metabolic 
fluxes. Similarly we should expect renewed attention to arise from the other term 
in our summation theorem, i.e. the one that refers to control by the diffusion. Of 
course the novelty lies not so much in the possibility of such control, but in the fact 
that it may account for a shortfall of metabolic control from the total of 1. 

We have also formulated the summation theorem for the case of inhomogeneity 
in both time and space, i.e. where a spatial gradient may be developing over time. 
In that case diffusion and time turn out to share in the total flux control of 1. 
In the case of homogeneity in time (i.e. at steady state), the former contribution 
disappears, in case of spatial homogeneity the latter. The experimental system 
that springs to mind is that of Belousov-Zhabotinsky reactions exhibiting time 
varying spirals (Zhabotinsky and Zaikin 1973), but closer to biology the wave type 
oscillations in yeast extracts (Mair et al. 2001). The most relevant application 
may be in developmental biology, but for this an extension to Hierarchical Control 
Analysis (Snoep et al. 2002) will be needed. 


While the hrst flux summation result, Eq. 11, and the hrst concentration sum¬ 
mation result, Eq. 16, are a simple extension of existing theory, the summation 
theorems concerning changes of spatial scale (Eqs. 18, 19, and 20) are fundamen¬ 
tally new, and introduce a concept that is new to Metabolic Control Analysis: 
control by size. 
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‘Control by the size’, the quantity measured by Cl, is a concept with some 
interesting properties. For instance, C£ and C£ can be both positive and negative. 
This is related to the fact that the size of a system can influence the ‘productive’ 
capacity of the system in different ways at the same time. Let us examine this for 
flux control: 

• Size can be an obstacle: if material has to travel over a distance L, as in the 
example of parallel membranes (Case 1 of Section 5), then an increase in L 
implies a decrease of the (average) concentration gradient, and therefore 
causes a reduction in diffusive flux. This is the dominant effect in the 
large-L limit of Case 1, as shown in Fig. 5, where it leads to a limit value 
oiCi of-1. 

• On the other hand, size can also be a resource: if a reaction takes place in 
a bulk region of linear size L, then L is also a measure of reactive capacity: 
an increase in size will result in an increase of total reactive flux. This is 
demonstrated by the second example, where the input- and output bound¬ 
aries are co-located and no material is forced to travel over a distance L. 
Here the control by size is positive, and in the limit of small L the control 
even approaches 1. 

The new summation theorems have a number of interesting implications. Per¬ 
haps the most striking one is that they impose an upper limit to the control that 
diffusion may have on the flux through the system. The first summation theorem 
suggested that this number be 1, i.e. if the biochemical and the transport processes 
were in excess and the system were large; the second, however, shows that it is only 
1/2 whenever the control by size is positive (i.e. when increasing the size results in a 
higer flux). Clearly, in such cases there must always be another process controlling 
flux, in addition to diffusion. In the example of Fig. 6 the other processes were the 
biochemical reactions. 

The example of Figure 2 also shows how changing the diffusion rate alters the 
concentration profiles. By increasing the diffusion rate we increase the width of 
a ‘boundary layer’, in which a significant part of the reactions takes place. As 
a result the flux is also increased. Note, however, that the relative increase in 
size of the boundary layer (which scales as A, since a; = ol'C = ^C) i® f^e square 
root of the increase of the diffusion coefficients (a^j = '^^)- This can be explained 
intuitively as follows. The size of the boundary layer is determined by counteracting 
forces; a thin boundary layer provides for fast spatial transport (on the scale of the 
layer thickness), but little reactive capacity, while for a thick layer the situation 
is inversed. An increase in diffusion rate is therefore used partially to increase 
transport rate, and partially to increase reaction capacity. The number 1/2 is the 
result of this trade-off. 
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